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Abstract 

In  this  paper,  we  develop  a  linearized  imaging  theory  that  combines  the 
spatial,  temporal  and  spectral  components  of  multiply  scattered  waves  as  they 
scatter  from  moving  objects.  In  particular,  we  consider  the  case  of  multiple 
fixed  sensors  transmitting  and  receiving  information  from  multiply  scattered 
waves.  We  use  a  priori  information  about  the  multipath  background.  We  use  a 
simple  model  for  multiple  scattering,  namely  scattering  from  a  fixed,  perfectly 
reflecting  (mirror)  plane.  We  base  our  image  reconstruction  and  velocity 
estimation  technique  on  a  modification  of  a  filtered  backprojection  method  that 
produces  a  phase- space  image.  We  plot  examples  of  point- spread  functions 
for  different  geometries  and  waveforms,  and  from  these  plots,  we  estimate  the 
resolution  in  space  and  velocity.  Through  this  analysis,  we  are  able  to  identify 
how  the  imaging  system  depends  on  parameters  such  as  bandwidth  and  number 
of  sensors.  We  ultimately  show  that  enhanced  phase-space  resolution  for  a 
distribution  of  moving  and  stationary  targets  in  a  multipath  environment  may 
be  achieved  using  multiple  sensors. 

(Some  figures  may  appear  in  colour  only  in  the  online  journal) 


1.  Introduction 

Traditional  wave-based  imaging  theory  makes  the  assumption  that  the  waves  travel  directly 
from  the  sensor  to  the  target  and  back.  However,  there  are  many  situations  in  practice  where 
complicated  scenes  cause  waves  to  multiply  scatter  from  non-target  objects  in  the  scene 
of  interest.  This  makes  imaging  targets  of  interest  more  difficult  because  of  the  additional 
unknowns  that  are  introduced.  A  number  of  researchers  have  considered  wave-based  imaging 
in  the  presence  of  multipathing  when  the  scene  of  interest  consists  of  stationary  targets 
illuminated  by  waves  transmitted  from  moving  platforms.  The  studies  [9]  and  [4]  showed  that 
a  backprojection  algorithm  that  exploits  multiply  scattered  fields  can  improve  image  fidelity 
if  it  is  possible  to  uniquely  identify  the  part  of  the  data  corresponding  to  each  scattering 
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wave  path.  In  particular,  this  work  showed  that  targets  oriented  in  certain  directions  can  be 
effectively  invisible  from  direct  scattering  but  may  become  visible  when  multipath  scattering 
is  used  in  the  image  formation  process.  The  studies  [2]  and  [5]  showed  that,  compared  to 
standard  backprojection,  a  backprojection  algorithm  that  exploits  multiply  scattered  fields  can 
improve  image  resolution  even  if  a  method  for  a  priori  separating  paths  in  the  data  is  not 
feasible.  This  paper  derives  a  data  model,  based  on  [1],  that  includes  the  effect  of  multiply 
scattered  waves  on  moving  targets  and  then  develops  a  corresponding  phase-space  imaging 
method  that  accounts  for  multiple  scattering  in  the  data.  We  consider  the  case  in  which  multiple 
sensors  interrogate  a  distribution  of  moving  targets  that  are  restricted  to  a  horizontal  plane. 
We  model  multiple  scattering  with  a  perfectly  reflecting  (mirror)  horizontal  plane. 

This  geometry  might  arise  in  a  variety  of  physical  problems,  including  ultrasound, 
microwave  tomography,  seismic  prospecting,  and  radar  and  sonar  imaging.  The  perfect 
reflecting  mirror  is  a  simple  model  for  the  ocean  surface,  for  example,  that  must  be 
considered  in  marine  seismic  surveys,  active  sonar  surveys  and  radar  imaging  in  marine 
environments.  Alternatively,  the  mirror  could  represent  a  reflecting  surface  in  ultrasound, 
microwave  tomography  or  radar  imaging.  In  the  discussion  below,  we  use  radar  terminology, 
but  the  theory  applies  equally  well  to  other  modalities. 

The  mathematical  model  is  discussed  in  section  2.  Imaging  is  addressed  in  section  3,  which 
outlines  the  backprojection-type  algorithm  for  forming  phase- space  images,  and  then  proceeds 
to  analyze  the  corresponding  point-spread  function  (PSF).  Section  4  outlines  simulations.  Many 
of  the  details  can  be  found  in  [6]  and  short  summaries  can  be  found  in  [7]  and  [8]. 


2.  Mathematical  data  model 

A  radar  receiver  collects  information  from  an  electromagnetic  wave  that  was  sent  from  a 
transmitting  antenna  and  that  subsequently  scattered  off  a  set  of  (stationary  or)  moving  targets 
in  a  scene.  In  this  work,  we  assume  that  the  scattered  waves  arrive  at  the  radar  receiver  via 
multiple  paths. 

We  assume  that  the  electromagnetic  field  propagates  in  the  lower  half-space  Q  = 
{(jci,  X2,  x^)  \x\,  X2  G  M,V3  <  z0}.  The  perfectly  reflecting  mirror  is  located  at  dQ  = 
{(xi,x2,x3)  \x\,x2  G  M,  V3  =  zo}.  We  denote  the  reflection  of  Q  with  respect  to  dQ 
by  Q*. 

Both  Q  and  Q*  consist  of  a  homogeneous,  isotropic  medium  in  which  the  wave  propagation 
speed  is  c.  The  geometry  is  shown  in  figure  1. 

We  make  the  assumption  that  one  component  x/r  of  the  electromagnetic  field  satisfies 

(V2  -  c_232  -  V(t,  X)  92)  ir(t,  x,y)  =  s(t  -  ty,  x,y),  (1) 

with  the  boundary  condition 

f(t,x,y)  |ieas2  =  0.  (2) 

In  order  to  define  V(t,x),  consider  first  an  object  located  at  position  x  that  moves  with 
velocity  v.  This  object  would  be  described  by  an  index-of-refraction  distribution  with  a 
parameter  v : 

Vv(t,i)  =  c-2[n2v(i-vt)-\\  (3) 

Multiple  objects  translating  with  the  same  velocity  v  can  also  be  described  in  the  form  (3); 
for  example,  multiple  cars  moving  along  a  highway  at  the  same  constant  velocity  v  could  be 
described  by  (3)  in  which  nv  consists  of  multiple  bumps,  each  bump  representing  one  car.  For 
objects  translating  with  velocity  v ,  we  write  qv(x)  =  c~2[n2{x)  —  1]. 
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Figure  1.  This  shows  the  geometry  with  one  transmitter  at  y  and  one  receiver  at  z.  The  mirror  is 
located  at  height  zq. 


Objects  moving  with  a  different  velocity  can  be  described  by  an  expression  of  the  form  (3) 
in  which  the  velocity  v  is  different.  To  describe  two-way  traffic  along  a  highway,  for  example, 
V  might  be  of  the  form 

V(t,x)  =  c~2\n2v(x-  vt)  -  l]  +c~2\n2_v{x+  vt )  -  l],  (4) 

where  nv  describes  the  cars  moving  in  one  direction  and  n_v  describes  the  cars  moving  with 
the  same  speed  in  the  opposite  direction. 

If  there  are  many  objects  moving  with  different  velocities,  then  we  can  describe  the 
distribution  of  multiple  moving  objects  as 

V(t,  x)  =  J  c~2q2{x  —  vt)  dv.  (5) 

We  see  that  qv  ( x )  can  be  considered  a  scattering  density  in  phase  space. 


2.1.  Incident  field 


The  incident  field  \j/m  is  the  field,  propagating  in  the  empty  reflecting  environment,  that 
emanates  from  a  transmitting  antenna  located  at  position  y.  The  field  \j/m  satisfies 


(V2  -  c~2dj)im  =Six-  y)s{t  -  ty ) 

(6) 

irm  ( t,x,y )  |an  =  0. 

(7) 

An  expression  for  the  incident  field  can  be  derived  by  using  the  method  of 
images.  Corresponding  to  the  source  location  y  e  we  define  the  virtual  source  at 
/  =  (yi,y2,2z0-y3). 

First,  the  outgoing  Green’s  function  satisfying 


(V2  -c-2dj)g{t,x,y)  =  S(x-y)m  (8) 

=  0  (9) 


is 


g(!,x,y)  =  gy(t,x,y )  -g?(t,x,y). 


where  gy  is  the  free-space  Green’s  function  with  a  source  aty: 


gy(t,x,y) 


*(*  -  M 

An\x-y\ 


(10) 

(11) 
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and  similarly  gy  is  the  free-space  Green’s  function  with  a  virtual  source  at  /.  Note  that  g 
vanishes  at  the  boundary  dQ. 

This  means  that  the  solution  to  (6)  with  (7)  is 


-  j git -T.DHi-yMT -<»)djdr 

4jt\x  —  y'\  4tt\x—  y\ 


(12) 


2.2.  Scattered  field 

We  split  the  field  into  the  incident  field  and  the  scattered  field  t/tsc: 

xlr  =  +  iAsc.  (13) 

We  subtract  (6)  from  (1)  and  use  (13)  to  arrive  at  the  partial  differential  equation  for  the 
scattered  field: 

(V2  -  c“2a2)  =  V (t,  x)d?ijfsc(t,  z,y) 

rc  lan  =  0.  (14) 


We  use  the  Born  or  single-scattering  approximation  [1],  replacing  \jr^c  on  the  right  side 
of  (14)  by  The  corresponding  solution  of  (14),  measured  at  the  receiver  y,  constitutes  the 
data  d.  With  the  Green  function  (10),  we  obtain 


d  (t,y,z)=  J 

-I 


4tv\z  —  x\ 


4jv\z  — 


4 7t  | z  —  x\ 


dv 


V(t\x)d;\lfin(t',x,  y)dt'dx 


+ 


-=fi~  J  <!•{*-  vt') 

r  8  (t r 

/  W 

rdt-t'-—)  r 


Sy  (t'  ty  1  ) 

v  /  A  +  '  AZ 


4n\x  —  y\ 


d  t  dx 


dv 


Sy  h  1  / 1 ) 

\  /  A S 


4jt\x  —/| 


dt  dx 


dv 


Sy  h  1  /') 

\  /  A-t-' 


dt  dx 


4tc\x  —  y  | 

(tf  -  ty  -  l-y—  ^ 

dv — ^ - - - —  dt  dx.  (15) 

4n\x-y\  v  7 


Equation  (15)  yields  a  model  for  the  scattered  wavefield.  Each  term  of  (15)  has  a  clear 
physical  interpretation.  The  first  term  (line  2)  corresponds  to  activating  the  transmitter  at  y 
with  a  waveform  sy,  starting  at  time  ty.  This  activation  results  in  a  wave  leaving  the  transmitter 
and  propagating  directly  to  the  target,  arriving  at  time  tr.  The  moving  target  that  at  time  t'  is 
located  at  x  was,  at  time  t  =  0,  located  at  x  —  vt' .  The  current  density  induced  on  the  target 
by  the  incident  field  is  assumed  to  be  proportional  to  the  second  derivative  of  the  transmitted 
signal;  this  is  a  result  of  the  Born  approximation.  After  interaction  with  the  target,  the  wave 
then  propagates  directly  from  x  to  z. 

The  second  term  of  (15)  (line  3)  corresponds  again  to  activating  the  transmitter  aty  with 
the  waveform  sy.  This  time,  however,  the  wave  reflects  from  the  mirror  on  its  way  to  the  target; 
by  the  method  of  images,  the  field  that  arrives  is  the  negative  of  the  same  wave  traveling  in 
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free  space  but  transmitted  from  the  virtual  location  /.  This  field  interacts  with  the  target  and 
propagates  to  the  receiver  as  before. 

The  third  term  of  (15)  (line  4)  corresponds  to  direct-path  propagation  between  the 
transmitter  and  the  target,  and  a  reflection  from  the  mirror  on  the  way  to  the  receiver; 
by  the  method  of  images,  the  received  wave  is  the  negative  of  the  wave  propagating  in 
free  space  that  would  have  been  received  at  the  virtual  receiver  location  z'.  The  fourth  term 
corresponds  to  mirror  reflections  between  the  target  and  both  the  transmitter  and  receiver. 

We  interchange  the  order  of  integration  and  perform  a  change  of  variables  x  =  x  —  vt'\ 


d  (t,y,z) 


/ 


s(t-t .  Sy(t'-ty- 


‘I 
‘I 
+  / 


4jt\x  +  vt'  —  z\ 

4jt\x  +  vt'  —  z\ 
8  (t  -  t’  - 


47T  \x  +  vt'  —  z 


J  qv  (x) 

\\ 

J  qv  (x) 

^  fqv 

-A\ 

— Sqv 


\x+vt'-y\ 


—  At'  dc  du 


4tx\x  +  vt'  —  y  | 

s, 

4tt\x  +  vt'  —y'\ 

i,  li±f=2l) 

4tt\x  +  vt'  —  y  | 

i, 


(x) 


(x) 


d tf  dx  dv 


dt'  dx  dv 


dtf  dxdv. 


(16) 


4tt\x  +  vt' —  z'\  J  4it  \x  +  vt'  —  y'\ 

Finally,  we  carry  out  the  t'  integration  in  (16)  by  using  properties  of  delta  functions  with  more 
complicated  arguments  as  in  [6].  In  particular,  in  the  first  term  of  (15)  (line  2),  for  example, 
we  find  that  the  delta  function  contributes  only  to  tr  satisfying 

\x  +  vt'  -  z\ 


0  =  t  -  t'  - 


(17) 


This  time  we  denote  tr  by  txz.  We  also  introduce  the  notation  Rxj(txz )  =  \x  +  vtxz  —  z\  and 
u  =  \f'z  (Ixz(t))  I  for  fX7  (V)  =  -1  -  Rxz(t’)/c. 

This  leads  to  the  general  data  set  for  a  single  source  y. 


d  (t,z,y) 


/ 


SyU^jt)  -ty-  (X) 

(4 tc  )  Rxz  (txz(t)^  RXy  (txz {t ))  {iXz,v(t) 


i 

i 

/ 


ty 


dxdv 


RXy'  yxz  (0  ) 


Mv  (x) 


(4lC  )  Rxz  (^jcz(O)  Rxy'  (^xz(O)  /4j <cz,v(t) 

Sy[ixAt)-ty-R-^^]qv  (X) 


(4tT)  Rxz '  (txz;  (/^))  RXy  (txzr  (^  ) )  fj^xz' ,  U  (0 


dxdv 


dxdv 


Sy  [tXzr  (t )  ty 


Rxy' fez' ( 0) 


Mv  (x) 


(4 n)2 Rxz' (ixAt))  Rxy'  {txz'  (0)  H'xz! ,v  (0 


dx  dv. 


(18) 


23.  Slow -moving  approximation 

We  consider  only  objects  that  move  significantly  slower  than  the  speed  of  light.  Consequently, 
we  assume  that  the  distance  traveled  by  the  objects  during  the  relevant  time  is  much  smaller 
than  the  distances  between  the  targets  and  the  transmitters  and  receivers.  Specifically,  we 
assume 

\v\t,  \v\2t2u>m!a/c  «;  \x-z\,\x-z!\,\x-y\,\x-y'\,  (19) 

where  comdx  is  the  effective  maximum  angular  frequency  of  the  signal  sy. 
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This  means  that  a  first-order  Taylor  approximation  may  be  made  such  that 

Rxzit)  =  \z-x+vt\=Rxz  (0)  +  Rxz  (0)  •  vt  +  O  ,  (20) 

where,  Rxz  (0)  =  x  -  z,  Rxz  (0)  =  \RXZ  (0)  |,  and  Rxz  (0)  = 

The  techniques  of  [3]  can  be  used  to  derive  a  time  dilation  or  a  Doppler  scale  factor 


1  - 


Rxy(0)-V 


1  + 


Rxzioyv 


(21) 


that  is  closely  related  to  the  Doppler  shift.  We  can  further  approximate  this  factor  as 


*xvyz 


=  1  —  (Kjrv(0)  +/?jrz(0))  • 


(22) 


We  immediately  see  that  the  Doppler  scale  factor  axvyx  depends  on  the  bisector  (Rxy(Q)  + 

Rxz(0)). 

With  the  approximations  (20)  and  notation  (21),  the  data  expression  (18)  for  the  slow- 
mover  case  becomes 


d s(z,y,t)  = 


/ 


View  •  (t  ~  ^  -  ty]qv  (x) 


(4 Tt)  RXZ  (0)  Rxy  (0  )  /ixx 


dxdv 


sy[axvyz  ■  (t  -  &f>)  -  M  _  ty\qv  {x) 

(47r )2 ~R~XZ  (0) ~R~dy  (0) [iXZ'V 

Sy[otXVyz'  ■  (f  -  )  -  ^7^  -  ty]qv  (x) 

(4tt)2  Rx/  (0)  Rxy  (0)  /w 

sy[axvy/  ■  (t  -  ^^)  _  -  ty]qv  (*) 

(4n)2  Rx-j  (0 )Rxy  (0)  w.i. 

In  order  to  more  succinctly  describe  the  data  model,  we  write 


/ 

/ 

/ 


dxdv 


dxdv 


dxdv. 


(23) 


Rxy  —  **10,  Rxy'  —  Rx  20,  Rxz  —  Rx  01?  Rxz'  —  Rx 02-  (24) 


In  other  words,  we  use  the  subscripts  v  jk  and  v  jk  as  placeholders  for  the  dependence  on  the 
object,  receiver  or  transmitter  location.  A  zero  is  added  in  the  respective  placeholder  if  there 
is  no  dependence  on  a  transmitter  and/or  receiver.  Moreover,  we  use  /,  k  =  1  for  the  actual 
sensors  and  j,k  =  2  for  the  virtual  sensors.  We  also  write  aVJ ^  as  avn,  etc. 

With  this  notation,  we  rewrite  equation  (23): 


where 


ds(z,y,0  =  f  X!  TA"  (*)  dxdu, 

J  j.k=  1 

(- 1  )j+k  s\axvjk  (t  -  M  _  ty  | 

(4n  )2  Rm  (0)  RXjo  (0) 


(25) 


(26) 


2.4.  Narrowband  waveform 

Most  radars  use  narrowband  waveforms,  which  may  be  written  as  s(t)  =  a(t)  e-1"0^,  where 
a(t)  is  slowly  varying  and  coojo  is  the  carrier  frequency  for  the  transmitter  y.  (Note  that  in 
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keeping  with  the  notation  defined  above,  the  zeros  indicate  that  there  is  no  object  or  receiver 
dependence.) 


J\k=  1  J 


exp  (-icooyo [axvjk  (t  -  ^ 


(' 


xa  \  t  —  ty  — 


RxOk  (0)  +  RxjO  (0) 


RxOk  (0)  RxjO  (0)  flxO k,v 

qv  (x)  dxdv. 


(27) 


where  we  use  the  approximation  avjk  =  1  in  the  slowly  varying  amplitude  a.  Substituting 
the  time-independent  exponential  phase  terms  into  the  function  (pXVjk  =  [Rxj0  (0)  + 
oivjkRxOk  (0)  +  cty],  we  rewrite  (27)  as 

d„<z,*o  =  £/(«)  V.y+‘ 


gi <Pxv  jk  g  lcox jO&v  jkt 


j,k=  1  * 


xa  \  t  —  tv 


RxjO  (0)  RxOk  (0)  l^x0k,v 

RxOk  (0)  -  Rxj0  (0)\  A  A 
1  qv  (x)  dx  dv 


=  J2  (~Oi+k  d 

j,k= 1 


(28) 


3.  Image  reconstruction 

In  this  section,  we  first  describe  how  we  form  an  image.  Then  we  describe  how  we  analyze  the 
reconstructed  image.  We  test  our  analysis  by  performing  numerical  experiments.  We  conclude 
with  a  summary  of  our  results. 


3.1.  Image  formation 


The  approach  we  take  to  form  an  image  of  our  simulated  data  involves  applying  a  matched 
filter  or  filtered  adjoint  to  the  data.  From  the  data  d (t,  z,  y),  we  form  an  image  I  ( p ,  u)  of  the 
objects  with  hypothesized  velocity  m,  that,  at  time  t  =  0,  are  located  at  position  p. 

Consider,  for  example,  the  case  where  the  transmitters  and  receivers  illuminate  only  one 
path,  say  the  direct  path.  In  this  case,  the  data  are  simply 


\t,z,y)  =  [Tnqv(x)](t,z,y)  =  J  Y Jj(x,v,t,z,y)qv(x) 


dxdv, 


(29) 


where  we  have  applied  (19),  (20)  and  the  narrowband  approximation  to  (26)  to  obtain 
Tv1  (x,  v,t,z,y ): 

(6>y\ 2  e^-e-iw  (  Rxz  (0)  Rxy  (0)  ^ 

1  —  I  - all  \  t-ty - - - £ -  )  .  (30) 


T"(*,  v,t,z,y )  = 


\4n/  Rxz  (0)  Rxz  (0)  /xA-it,  \  c 

In  order  to  form  an  image  I(p,u )  at  position  p  and  velocity  u ,  we  apply  the  adjoint 
operator  to  d^1  (t,  z,  y): 


up,  u)  =  vn*  41 

=  X!  f  tn*(P’  u,t,z,y)d„(t,z,y)dt 

Z,y  ^ 


T ju(p,  u,t,z,y)  y“(jc,  v,  t,z,y)Atqv{x)Axdv 
( p ,  x,  u,  v)qv(x)  dxdu, 


(31) 
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where  the  kernel 

IC, 


"O,  X,  u,  v)  =  f  u,  t,z,y)Tkl(x,  v,t,z,y)dt 

J  z,y 


(32) 


is  the  PSF  for  the  imaging  system.  We  note  that  the  PSF  and  the  corresponding  image  are 
in  general  complex- valued.  The  Cauchy-Schwartz  inequality  implies  that  the  integral  kernel 
(32)  is  maximized  when  (p,  u)  =  (x,v).  In  particular,  each  term  is  bounded  in  magnitude  by 


I  / 


Y x(p,u,t,  v,t,  .  ..)d  t 


<  IIT/vO*’  u,  . .  - )  1 1 2 1 1  (  a:  ,  v,  . .  .)||2, 


(33) 


where  ||Y(. . .) II2  =  (/  | Y(. . .  ,t,  . .  .)|2  and  equality  is  attained  in  (33)  when 

Y^1  ( p ,  u,  t,  z,  y)  =  Y^1  (x,  v,t,z,y).  This  occurs  when  ( p ,  u)  =  (x,  v). 

If  we  can  identify  the  parts  of  the  data  that  correspond  to  the  different  paths,  then  we 
can  form  separate  images  for  each  path  and  coherently  (or  noncoherently)  add  the  resulting 
images: 

/(/>,  u)  =  E  Vij*diJ.  (34) 


hj 

If  it  is  not  possible  to  distinguish  paths,  then  the  data  are  of  the  form  d  =  ^  dlm ,  we  apply 
the  adjoint  operator  =  J2i  j  'Pn*  t0  form  an  image: 

I(p,  u)  =  VydN(y,z,t) 

y  Y x(p,  u,  t ,  z,  jOT/v(x,  ^  z,  y)  dAyv(x)  dxdv,  (35) 

where  the  subscript  N  denotes  the  narrowband  slowly  moving  case.  The  image  that  is  formed 
gives  us  an  approximation  to  the  true  reflectivity  function  qv  (jk;),  added  to  copies  of  qv  in  the 
wrong  location. 

3.2.  Image  analysis 

In  order  to  analyze  the  image  in  the  case  when  data  from  the  different  paths  cannot  be  separated, 
we  need  to  investigate  the  relationship  between  the  image  and  the  true  reflectivity  qv  (. x )  . 

To  do  this,  we  use  (28)  in  (35).  The  result  is 

2 

I(p,u)=  y  (— 1  )j+kKJNk(p,  u;  x,  v)qv(x)  dxdv,  (36) 

j,k=  1 

where  the  sum  of  all  terms  is  the  PSF  for  this  imaging  system.  Each  is  given  by 

(Rpok(0)  +  RPjo(0))' 


kn  {p-  X,u,v)  =  -  j  a>lj0a*k  ^ 

V 


X  a  jk  I  t  ty 


(RjcOk  (0  )  +  ^/0(0)) 


xe 


~1(Pp jk  p'^0  17O Olpjk  pl Vx jk  p -kOJQ jOO(xjk 


d  t. 


R  pOk  (0  )R  p  j0  /X  ^0^,  u 
RxOk  (0  )Rx jOF'xOk,  u 

Equation  (37)  can  be  further  simplified  by  noting 

{ Pxvjk  (Ppujk  =  ^  \R-xj0  (0)  ^v  jkRxOk  (0)  T"  Cty  (RpjQ  (0)  &ujkRpQk  (0)  T"  Cty )  ] 

=  — ^[^jcyo(O)  “  RpjoiP)  —  OtvjkRxOk(0)  +  <XujkRp0k(®)] 

=  ^[RxjoiO)  -  Rpjo(0)  -  (1  +  pvjk)RX0k(Q)  +  (1  +  PujkWM 0)] 

=  -^p-[cA*xpjk  —  kv  jkR*0k  (0)  +  PujkRpOk(0)], 


(37) 


(38) 
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where 


^‘T'xpjk  — 


Rxjo(Q)  -  Rpjo(0)  -  Rxok(Q)  +RPok(Q) 


(39) 


With  this  notation,  we  can  write  (37)  as 


4k(p 


f  2  /  fi 

,  U.  X,  V  )  —  /  COqi^Aj^  (cOqjq [fipjk  fixjk]  •>  ^TXpjk)&k P  y^Ojol^pjk  fxjk\ 


Rxjo(O) 


xexp  [  1^M[Rxj0(0)  -  Rpj0(0)] 


RpOk(0)RpjO  (0) l^p0k,u 
Rx0k(^)Rxj0  (®)  Tx0k,u 


d  t, 


(40) 


where  avyz  ~  1  +  fivyz  (see  equation  (22))  and 


A(S),  t)  =  e_ko°r  J  a*  (t  -  r )a(t)  e3*  d t  (41) 

is  the  well-known  narrowband  ambiguity  function  [10]. 

We  recall  that  the  ambiguity  function  characterizes  the  ability  of  a  given  waveform  to 
provide  range  and  (down-range)  velocity  information.  Thus,  formulas  (36)  and  (40)  provide 
a  way  to  examine  the  performance  of  a  system  for  imaging  moving  targets  in  terms  of  the 
sensor  locations  and  the  waveforms  they  transmit.  These  formulas  thus  provide  a  theoretical 
basis  for  the  allocation  of  radar  resources. 


4.  Numerical  examples 

In  this  section,  we  provide  some  numerical  examples  that  show  how  the  number  of  sensors 
affects  the  image  reconstruction  of  a  stationary  and  a  moving  target.  We  analyze  the  resolution 
of  the  data  by  plotting  a  PSF,  which  is  equivalent  to  an  image  of  a  single  point-like  target.  All 
simulations  used  the  following  parameters. 

(i)  Target.  The  target  is  a  single  point-like  scatterer  located  in  the  horizontal  plane.  For  the 
stationary  target  images,  the  true  target  is  located  at  x  =  (25  m,  25  m)  in  a  scene  in  which 
the  v  axis  runs  from  0  to  50  m  and  the  y  axis  runs  from  0  to  50  m.  For  the  moving  target 
images,  the  target  is  located  at  (35  m,  35  m)  in  a  scene  where  both  v  and  y  axes  run  from 
0  to  70  m. 

(ii)  Waveform.  Each  sensor  transmitted  a  known  stepped-frequency  waveform  with  218  pulses 
stepping  through  the  frequencies  of  30-40  GHz  uniformly.  Alternatively,  this  can  be 
thought  of  as  a  signal  with  a  35  GHz  center  frequency  and  a  10  GHz  bandwidth. 

(iii)  Sensor  geometry.  The  sensors  lie  around  the  scene  of  interest.  See  figure  2  for  the  exact 
placement.  Note:  the  velocity  plots  use  only  the  nine  transmitters  and  one  receiver  set-up. 

(iv)  Reflector.  The  reflector  was  a  perfectly  reflecting  mirror  situated  100  m  above  the  scene 
of  interest.  The  location  of  the  reflector  is  assumed  to  be  known. 

In  all  cases,  the  images  are  spatial  images.  The  units  on  both  the  horizontal  and  vertical 
axes  are  in  meters. 

Figures  3  and  4  show  reconstructions  of  a  stationary  target;  figure  3  shows  the  case  in 
which  we  are  able  to  identify  the  part  of  the  received  signal  coming  from  different  transmitters, 
whereas  figure  4  shows  the  case  in  which  we  do  not  have  this  information.  All  the  images  in 
figures  3  and  4  show  the  spatial  image  at  velocity  (0,  0). 

We  note  that  since  figure  3  is  reconstructed  using  prior  information  about  the  location 
of  the  reflecting  mirror,  and  consequently  the  wave  path  is  known,  there  are  no  multipath 
ambiguities  because  the  algorithm  matches  the  correct  time  delay  with  the  correct  path. 

As  we  increase  the  number  of  sensors,  the  resolution  improves:  the  red  region 
corresponding  to  the  target  is  thinner  than  the  backprojection  ellipses.  Since  each  image 
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Figure  2.  Geometry  of  the  sensor  placement.  Top  image:  three  transmitters  and  one  receiver. 
Bottom  image:  nine  transmitters  and  one  receiver. 


cell  represents  2  m  of  physical  distance,  we  estimate  the  resolution  for  the  stationary  point 
target  in  this  example  to  be  about  2  m. 

In  figure  4,  on  the  other  hand,  the  image  is  reconstructed  without  the  knowledge  of  the 
paths  traveled  by  different  parts  of  the  received  signal.  As  the  figure  shows,  multipath  artifacts 
now  appear.  We  note  that  there  is  a  relationship  between  the  configuration  of  the  sensor 
positions  and  the  artifacts.  Adding  more  sensors  diminishes  the  artifacts. 

In  figure  5,  the  target  was  moving  with  velocity  (3,  3)  x  10  m  s-1.  Figure  5  shows  spatial 
images  at  two  different  velocity  hypotheses,  the  velocity  (1 ,  3)  x  10  m  s-1  (top)  and  the  correct 
velocity  (bottom).  The  numerical  examples  presented  in  this  paper  are  only  a  representative 
subset  of  velocity  images.  The  reader  is  referred  to  [6]  for  a  larger  set  of  velocity  images.  We 
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Figure  3.  Stationary  object  PSF — perfect  reflector,  known  paths.  Top  image:  three  transmitters 
and  one  receiver.  Bottom  image:  nine  transmitters  and  one  receiver.  The  vertical  scale  (image 
amplitude)  is  the  same  for  both  plots. 


can  determine  the  correct  velocity  by  choosing  the  image  with  the  best  target  focus.  In  this 
way,  we  may  use  our  modified  backprojection  method  to  estimate  the  velocity  of  a  moving 
target.  In  figure  5,  we  show  that  if  we  use  a  priori  information  about  the  geolocation  of  a 
target,  then  we  may  estimate  that  target’s  velocity  by  examining  the  four-dimensional  PSF. 
As  we  see  in  figure  5,  the  target  is  geolocated  with  the  best  resolution  at  the  correct  velocity 
(3, 3)  x  10  ms-1.  At  the  incorrect  velocity,  the  true  target  disappears  and  artifacts  appear  at 
the  wrong  locations. 

From  our  velocity  sampling,  we  estimate  the  resolution  in  velocity  to  be  about 
10  m  s-1,  and  the  image  shows  the  spatial  resolution  to  be  about  3  m.  In  this  case,  the 
waveform  we  used  provided  better  spatial  resolution  than  velocity  resolution.  A  waveform 
with  better  Doppler  resolution  would  provide  better  velocity  resolution  but  presumably  poorer 
spatial  resolution.  There  is  currently  no  theory  that  will  explicitly  quantify  the  relationship 
between  the  waveform,  the  sensor  locations  and  the  phase-space  resolution.  This  is  left  as  an 
open  problem. 
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Figure  4.  Stationary  object  PSF — perfect  reflector,  unknown  paths.  Top:  three  transmitters  and  one 
receiver.  Bottom:  nine  transmitters  and  one  receiver. 
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Figure  5.  Spatial  images  of  a  moving  object  for  different  velocity  hypotheses.  This  image  is 
zoomed  to  show  the  spatial  range  [20  m,  60  m]  x  [20  m,  60  m].  Top:  velocity  (1,3)  x  10  ms-1 
(incorrect).  Bottom:  velocity  (3,  3)  x  10ms_1  (correct). 
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5.  Conclusions 

In  this  paper,  we  develop  a  method  for  forming  a  phase- space  image  of  multiple  moving 
objects  in  a  known  multipath  environment  from  scattering  data  received  at  multiple  sensors 
from  multiple  transmitters  that  could  be  transmitting  different  waveforms.  The  image  is  formed 
via  a  modified  backprojection  algorithm.  In  particular,  the  approach  we  take  to  arrive  at  the 
image  formation  algorithm  involves  developing  a  physics-based  data  model  for  the  information 
we  expect  to  receive  at  the  receiving  sensor  and  subsequently  applying  an  adjoint  to  the  data 
model  in  order  to  form  the  image.  The  resultant  phase-space  image  is,  in  general,  a  six¬ 
dimensional  image;  in  the  case  when  objects  and  their  motion  are  restricted  to  a  plane,  the 
phase- space  image  is  four  dimensional. 

We  characterize  the  performance  of  the  imaging  system  by  its  PSF,  for  which  we  obtain  an 
explicit  formula.  Numerical  exploration  of  the  PSF  provides  information  about  the  degree  to 
which  the  phase-space  image  contains  information  about  position  and  velocity.  In  particular, 
from  plots  of  the  PSF,  we  are  able  to  estimate  the  four-  or  six-dimensional  resolution  for  the 
image  and  ascertain  that  multiple  sensors  do  indeed  improve  the  image  resolution. 

These  results  show  that  the  modified  backprojection  algorithm  developed  in  this  paper  is 
able  to  recover  position  and  velocity  information  from  moving-object  scattering  data  collected 
in  a  known  multiple- scattering  environment.  We  also  find  that  the  PSF  provides  an  insight 
into  system  design  specifications,  such  as  how  many  sensors  are  needed  and  what  waveforms 
should  be  transmitted  in  order  to  produce  a  certain  resolution  image. 

Much  remains  to  be  done.  Our  approach  has  assumed  a  simple,  known  multipath 
environment,  perfect  clock  synchronization  for  coherent  imaging  and,  in  many  cases,  transmit 
waveforms  that  can  be  distinguished  from  one  another  by  the  receivers.  We  have  also 
used  a  simple  point-like  target  scattering  model.  Other  imaging  approaches,  such  as  those 
incorporating  sparsity  and  micro-Doppler  techniques,  remain  to  be  explored. 
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